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ABSTRACT 

The  numerical  prediction  of  the  aerodynamic  noise  radiated 
by  an  isolated  airfoil  is  performed  using  a Computational 
AeroAcoustics  (CAA)  method.  This  hybrid  method  combines  (i) 
a simulation  of  the  nearfield  unsteady  flow  and  (ii)  an  acoustic 
method  to  estimate  the  noise  radiated  in  the  far  field.  This 
process  is  applied  to  a symmetrical  NACA0012  airfoil  with  a 
constant  section  and  a blunted  trailing  edge  (TE),  at  a Mach 
number  of  0.205  and  an  angle-of-attack  of  5°.  The  Reynolds 
number  based  on  the  airfoil  chord  is  2.86  millions.  The 
computational  domain  has  a spanwise  extent  representing  3.3  % 
of  the  chord. 

The  nearfield  unsteady  flow  around  the  airfoil  is  computed 
via  a compressible  three-dimensional  Large  Eddy  Simulation 
(LES).  This  flow  exhibits  a superimposed  pressure  fluctuation 
field  which  presents  the  qualitative  and  quantitative  features  of 
the  TE  noise  generated  by  the  acoustic  scattering  of  (i)  the 
turbulent  boundary  layers  converted  on  both  airfoil  sides 
(broad-band  noise)  and  (ii)  the  alternated  vortex  shedding 
generated  by  the  TE  bluntness  (narrow  band  component).  Due  to 
the  strong  stretching  of  the  LES  computational  grid,  which  acts 
as  an  acoustic  low-pass  frequency  filter,  this  acoustic  field 
cannot  radiate  farther  than  a half-chord  from  the  body. 

Consequently,  the  LES  must  be  relayed  by  an  acoustic 
propagation  method  to  correctly  simulate  the  farfield  noise.  The 
most  readily  available  methods  are  integral  techniques  which 
only  need  data  on  a control  surface  or  in  a control  volume.  A 
former  estimation  is  performed  via  the  Kirchhoff  integration  of 
pressure  data  on  a control  surface  enclosing  the  airfoil  and  its 
wake.  The  position  of  this  surface  is  subjected  to  a careftil 
parametric  study.  Since  very  large  spanwise  integral  lengh 
scales  are  observed  on  this  surface,  the  acoustic  field  is  assumed 
two-dimensional  there,  leading  to  interesting  qualitative  results 
but  overestimated  noise  levels.  On  the  other  hand,  wall  pressure 
spanwise  integral  lengh  scales  reasonably  match  experimental 
values,  suggesting  to  use  the  Ffowcs  Williams-Hawkings 
integration  of  wall  pressure  data.  This  method  provides  farfield 
noise  levels  which  agrees  well  with  experimental  data. 

Both  methods  assume  that  the  mean  flow  is  uniform  outside 
the  control  surface,  which  slightly  distorts  the  solution  since  the 
actual  flow  is  non-homogeneous.  In  order  to  avoid  such 
distorsion,  one  can  resort  to  a numerical  method  based  on 
discretized  linearized  Euler  equations  (LEE)  in  which  the  total 
field  is  splitted  into  a mean  flow  and  a pertubation  field.  The 
E3P  code  (Propagation  of  small  perturbations  using  Euler 
equations)  is  currently  developed  for  that  purpose.  The  present 


paper  shows  its  application  to  several  complex  acoustic 
scattering  problems,  then  it  shows  how  this  code  could  be 
coupled  to  LES  calculation  on  a two-dimensional  interface. 

NOTATIONS 

ao  : sound  celerity  in  the  fluid  at  rest, 

A : characteristic  length  scale, 

Fj  : inviscid  fluxes, 

Fu  : viscous  fluxes, 

fea  : selective  function  fUl  (0) , 

® : wavenumber-frequency  spectrum 

H : trailing  edge  thickness, 

k : wavenumber, 

K : local  instantaneous  fluid  thermal  conductivity, 

*0  : thermal  conductivity  of  the  fluid  at  rest, 

K '■  small  scale  kinetic  energy, 

\ '■  Acoustic  wavelength  Ao  = coSc0, 

K ■ flow  Mach  number  Af0  = U0la0 , AtJ  = lf/yRT0 , 

A1  : local  instantaneous  dynamic  fluid  viscosity, 

Po  : dynamic  viscosity  of  the  fluid  at  rest, 

A : turbulent  dynamic  viscosity, 
n : unit  vector  normal  to  5,  oriented  towards  outside, 

© • frequency, 

p '■  local  instantaneous  pressure  in  the  flow. 

Pa  : pressure  of  the  fluid  at  rest, 

p'  ■ pressure  fluctuation,  identified  with  acoustic  pressure 
in  the  far  field  p’  = p-  pQ, 

Qc  '■  conservative  flow  variable, 

P : local  instantaneous  density  in  the  flow, 

Pa  '■  density  of  the  fluid  at  rest, 

P,  : density  fluctuation  in  the  flow  p'=  p-p0, 

Pr  : Prandtl  number  Pr  = p0Cp  / Ka 
A : Reynolds  number  A = PaUaC, I po  , 

Sj  : filtered  shear  stress  tensor 

S '■  Kirchhoff  integration  surface, 

TE  : trailing  edge. 
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TBL  : turbulent  boundary  layer, 

X '■  observer  location  in  the  airfoil  framework, 

Y ■ source  location  in  the  airfoil  framework, 
t : reception  time, 

T : emission  time, 

T '■  local  instantaneous  fluid  temperature, 

T0  : temperature  of  the  fluid  at  rest, 
u : local  instantaneous  fluid  velocity, 

U0  : upstream  flow  velocity, 

: local  mean  flow  velocity  outside  the  boundary  layer. 


INTRODUCTION 

The  aerodynamic  noise  generated  by  the  high  lift  devices  - 
HLD,  slats  and  flaps  - of  large  airliners  is  an  important 
contributor  to  the  total  radiated  airframe  noise,  especially  in 
approach  configuration.  It  is  widely  admitted  that  this  noise 
source  will  have  to  be  reduced  for  the  next  generation  of  aircraft 
to  meet  the  incoming  strengthened  noise  regulations  in  the 
vicinity  of  airports.  This  goal  will  be  reached  by  the  design  of 
new  HLD  concepts  incorporating  specific  noise  reduction 
devices.  Although  experiment  will  still  play  an  important  role  in 
the  design  of  such  devices,  it  is  thought  that  significant  results 
will  be  obtained,  with  lower  costs  and  shorter  delays,  from  the 
numerical  simulation,  particularly  considering  the  spectacular 
continuing  progress  of  Computational  AeroAcoustics  (CAA) 
applications. 

At  the  present  time,  the  problem  of  the  numerical  simulation 
of  HLD  noise  is  still  beyond  the  capabilities  of  complete  Direct 
Numerical  Simulation  (DNS).  It  is  most  likely  that,  in  the  next 
future,  only  hybrid  methods,  in  which  the  nearfield  turbulent 
flow  and  the  farfield  noise  are  computed  separately,  will  be  able 
to  solve  complex  aeroacoustic  problems  involving  HLD. 

Several  Computational  Fluid  Dynamics  (CFD)  techniques 
are  used  in  CAA  hybrid  methods.  The  less  computationally 
expensive  methods  are  based  on  steady  Reynolds-Averaged 
Navier-Stokes  (RANS)  computations,  in  conjonction  with 
stochastic  models  of  the  wavenumber-frequency  spectrum  of  the 
turbulence  '.  RANS  methods  are  also  used  to  compute  unsteady 
flows  around  airfoils  for  aeroacoustic  predictions  IJ,  but  Large 
Eddy  Simulation  (LES)  is  often  considered  as  one  of  the  most 
promising  method  in  the  context  of  CAA  applied  to  HLD 
noise 

On  the  other  hand,  a large  panel  of  acoustic  techniques,  with 
various  complexity  and  capabilities,  can  be  used  for  the 
prediction  of  farfield  noise.  The  most  practical  formulations 
belong  to  the  family  of  integral  methods  such  as  (i)  acoustic 
analogies,  including  Lighthill's  analogy  5,7  and  the  Ffowcs 
Williams-Hawkings  (FW-H)  equation  3J'\  (ii)  the  Boundary 
Element  Method  (BEM)  9 and  (iii)  the  Kirchhoff  integration 
method. 

Most  integral  methods  assume  that,  beyond  a given  distance 
from  the  noise  sources  and  body  surfaces,  the  sound  propagates 
in  a medium  at  rest,  or  moving  with  uniform  velocity.  This 
assumption  may  become  a strong  limitation,  especially  when  the 
radiated  noise  results  from  a surface  integration  on  a control 
surface  which  is  located  near  solid  walls  where  velocity 
gradients  are  significant.  In  that  case,  only  the  discretized  LEE, 
the  equations  governing  the  acoustic  propagation,  may  account 
for  the  propagation  in  non-homogeneous  flows.  Although  Euler 


equations  have  been  widely  used  in  CFD  applications  in  the 
past,  this  know-how  is  of  little  help  for  acoustic  simulation, 
since  higher  order  schemes  are  needed  to  ensure  numerical 
accuracy  and  low  dispersion  of  acoustic  waves  l0.  At  the  present 
time  hybrid  CAA  methods  combining  unsteady  CFD  and 
discretized  LEE  have  only  been  applied  to  free  field  flows  like 
jets  and  mixing  layers,  for  which  cartesian  grids  make 
computations  particularly  convenient.  Their  application  to 
realistic  geometries  including  airfoils  needs  curvilinear  grids 
on  which  the  use  of  high  order  finite  difference  schemes  is  not 
straightforward. 

First  CAA  studies  at  ONERA  focused  on  trailing  edge  (TE) 
noise  of  isolated  airfoils,  a basic  aeroacoustic  mechanism  which 
is  present  in  complex  HLD  geometries.  In  an  former  study  ', 
only  the  uncompressible  unsteady  flow  at  the  vicinity  of  the  TE 
was  simulated  via  LES,  which  was  valid  since  aeroacoustic 
sources  of  such  an  airfoil  are  known  to  be  concentrated  in  this 
region.  For  this  purpose,  the  TE  was  modelled  by  a thick  flat 
plate  with  an  imposed  mean  velocity  profile  of  a realistic 
boundary  layer  on  each  side.  This  simplification  was  also 
governed  by  the  capabilities  of  PEGASE,  the  used  LES  code 
which  only  handled  cartesian  meshes.  Then  the  radiated  noise 
was  obtained  using  the  formulation  Ffowcs  Williams-Hall  11 
based  on  Lighthill’s  acoustic  analogy  and  the  acoustic  Green's 
function  of  an  infinite  half-plane  (an  approach  that  has  been 
recently  followed  by  Wang  6).  Despite  encouraging  results,  this 
exploratory  study  presented  two  major  weaknesses  since  (i)  the 
exact  airfoil  shape  was  not  correctly  taken  into  account  in  the 
mechanism  of  acoustic  diffraction  of  quadrupolar  sources  and, 
(ii)  the  flow  simulation  did  not  comprise  the  whole  airfoil. 

In  a more  recent  study  \ point  (i)  was  addressed  by  the  use 
of  an  innovative  acoustic  method,  also  based  on  Lighthill's 
acoustic  analogy,  but  associated  to  a Boundary  Element  Method 
to  compute  the  exact  acoustic  Green’s  function  of  the  studied 
airfoil.  The  method  was  restricted  to  cylindrical  infinite  surfaces 
with  constant  sections,  for  which  the  determination  of  the  3D 
acoustic  Green's  function  is  narrowed  to  a much  simpler  2D 
problem  and  to  acoustic  radiation  in  far  field.  The  application  of 
this  method  to  more  realistic  airfoil  shapes  (a  NACA  0012 
symmetrical  profile  and  a high  lift  wing  section  with  slat  and 
flap),  but  with  non-realistic  unsteady  flow's,  showed  that,  as  long 
as  the  unsteady  three-dimensional  flow  is  known  in  the  vicinity 
of  a scattering  surface  of  any  shape,  the  resulting  radiated  noise 
can  be  predicted.  Consequently,  this  method  is  particularly 
attractive  when  the  turbulent  flow  is  computed,  for  example  by 
using  LES,  in  a domain  which  covers  only  a small  fraction  of 
the  whole  rigid  surface. 

The  present  paper  addresses  point  (ii),  the  objective  being  to 
perform  an  acoustic  prediction  based  on  the  LES  of  the  unsteady 
flow  around  a 3D  lifting  airfoil  with  a blunted  TE. 

The  paper  is  organized  as  follows.  The  first  part  of  the  paper 
is  devoted  to  the  unsteady  flow  simulation  via  LES.  The 
theoretical  background  of  LES  and  its  implementation  in  an  in- 
house  industrial  CFD  software  are  described.  Then  the  airfoil 
geometry  and  the  computational  grid  are  presented,  as  well  as 
the  numerical  procedure  and  storage  method.  A few  basic 
aerodynamic  results  are  given  at  the  end  of  this  part. 

In  the  second  part,  the  pressure  field  computed  via  LES  near 
the  airfoil  is  analyzed  in  details.  It  exhibits  the  qualitative  and 
quantitative  features  of  an  acoustic  field  noise,  generated  by  the 
acoustic  scattering  at  the  airfoil  TE  of  (i)  the  turbulent  boundary 
layers  convected  on  both  sides  (broad-band  noise)  and  (ii)  the 
alternated  vortex  shedding  generated  by  the  TE  bluntness 
(narrow  band  component).  It  is  shown  that,  due  to  the  strong 
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stretching  of  the  LES  computational  grid,  which  acts  as  an 
acoustic  low-pass  frequency  filter,  this  acoustic  field  cannot 
radiate  farther  than  a half-chord  from  the  body.  Consequently, 
the  LES  must  be  relayed  by  an  acoustic  propagation  method  to 
correctly  simulate  the  farfield  noise.  The  most  readily  available 
methods  are  integral  techniques  which  only  need  data  on  a 
control  surface  or  in  a control  volume. 

The  third  part  of  the  paper  presents  a noise  estimation 
performed  via  the  Kirchhoff  integration  of  pressure  data  on  a 
control  surface  enclosing  the  airfoil  and  its  wake.  The  position 
of  this  surface  is  subjected  to  a careful  parametric  study.  Since 
veiy  large  spanwise  integral  lengh  scales  are  found  on  this 
surface,  the  acoustic  field  is  assumed  two-dimensional  there, 
leading  to  interesting  qualitative  results  but  overestimated 
farfield  noise  levels. 

On  the  other  hand,  wall  pressure  spanwise  integral  lengh 
scales  reasonably  match  experimental  values,  which  suggested 
to  use  a Ffowcs  Williams-Hawkings  integration  of  wall  pressure 
data.  This  application  is  presented  in  the  fourth  part  of  the  paper. 
This  method  provides  farfield  noise  levels  which  agrees  well 
with  experimental  data. 

However,  both  integral  methods  assume  that  the  mean  flow 
is  uniform  outside  the  control  surface,  which  may  distort  the 
solution  if  the  actual  flow  is  non-homogeneous.  In  order  to 
avoid  such  distorsion,  one  has  to  resort  to  a method  based  on 
discretized  LEE  in  which  the  total  field  is  splitted  into  a mean 
flow  and  a pertubation  field.  The  E3P  code  (Propagation  of 
small  perturbations  using  Euler  equations)  is  currently 
developed  for  that  purpose.  The  fifth  part  of  the  paper  shows  its 
application  to  several  complex  acoustic  scattering  problems, 
then  it  shows  how  this  code  could  be  coupled  to  LES  calculation 
on  a two-dimensional  interface. 


LARGE  EDDY  SIMULATION 


Code  implementation 

The  LES  method  described  below  has  been  implemented  in 
the  FLU3M  solver,  an  industrial  CFD  software  which  has  been 
developed  at  ONERA  for  several  years.  This  solver  is  based  on 
the  discretization  of  the  compressible  Navier-Stokes  equations 
on  multi-block  structured  meshes  by  a finite  volume  technique. 
A second  order  accurate  implicit  temporal  integration  is 
achieved  thanks  to  an  approximate  Newton  method.  General 
informations  about  this  solver  are  available  in  refs.  ” “ 


Governing  equations 
Filtered  Navier-Stokes  equations 

A dimensionless  form  of  the  Navier-Stokes  equations  is 
retained,  which  means  that,  in  equations  (1)  to  (8),  the  local 
instantaneous  velocity  iT,  temperature  T,  density  p,  abscissa  x , 
dynamic  viscosity  p and  thermal  conductivity  K are  normalised 
respectively  by  U0 , T0 , p0 , C,  p^  and  ic0 . 

The  three-dimensional  unsteady  filtered  Navier-Stokes 
equations  are  used  for  a viscous  compressible  Newtonian  fluid. 
Any  flow  variable  <p  can  be  written  as  <p-  p +<p' , where  <j> 
represents  the  large  scale  part  of  the  variable  and  <p'  its  small 


scale  part.  The  filtering  operator,  classicaly  defined  as  a 
convolution  product  on  the  computational  domain  IS,  is  assumed 
to  commute  with  time  and  spatial  derivatives.  Moreover  it  is 
convenient  for  the  clarity  of  the  equations  to  introduce  the  Favre 
filtering  <j>  = ptpj p . In  conservative  form,  the  filtered  Navier- 
Stokes  equations  can  be  expressed  in  three-dimensional 
Cartesian  coordinates  (x,  ,x2,x3)  as  : 


dQ<  , dFi 

dr  dx. 
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where  Qc , Ft  and  (assuming  some  approximations  17  and  using 
Boussinesq  eddy  viscosity  hypothesis)  Fj  are  defined  by  : 
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where  the  eddy  viscosity  p,  must  be  expressed  by  a subgrid 
scale  model.  These  equations  are  supplemented  with  the  filtered 
equation  of  state : 


Suberid  scale  modelling 
Selective  mixed  scale  model 


(6) 


The  selective  mixed  scale  model,  developed  by  Sagaut  and 
Lenormand,  has  been  retained,  because  it  realizes  a good 
compromise  between  accuracy,  stability  and  computational  cost 
l7'".  More  particularly,  the  use  of  a selective  function  allows  to 
handle  transitional  flows  '9,  which  is  one  of  the  keypoint  of  the 
present  application. 


The  eddy  viscosity  is  given  by  the  non-linear  combination  : 

p,  = p/8„  (e^Ajo.ss)/^  (r<)/Z  (7) 

of  the  filtered  shear  stress  tensor  Sy(ft),  a characteristic  length 
scale  A (given  by  the  cell  volume  of  the  mesh),  the  small  scale 
kinetic  energy  Kc  = t/o.5  («,  )(«,  ) , a constant  parameter  C„  = 
0.06  (derived  in  reference  "’)  and  the  selective  function  : 

/a,(0)  = 


f 1 if  0 >60  = 10° 

[ tan4  (6/2  y tan4  (9uf  2)  otherwise 


(8) 


where  9 is  the  angle  between  instantaneous  and  filtered 
vorticity.  The  test  field  (tt, ) is  extracted  from  the  resolved  field 
(«,•)  = ii,  -Uj  , employing  an  averaging  test  filter  denoted  by  : 
5/  = 0.25u(_|  + 0.5f/,  + 0.25u  i+1 . 
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Numerical  method 

The  Navier-Stokes  equations  are  discretized  using  a cell- 
centered  finite  volume  technique  and  structured  multi-block 
meshes.  The  viscous  fluxes  are  discretized  by  a second-order 
accurate  centered  scheme.  For  efficiency  reason,  an  implicit 
time  integration  is  employed  to  deal  with  the  very  small  grid 
size  encountered  near  the  wall.  Indeed  the  large  disparity 
between  the  acoustic  wave  speed  and  the  advection  velocity  at 
low  Mach  number  renders  explicit  time  integration  unefTicient. 
This  is  due  to  the  fact  that  numerical  stability  considerations 
impose  small  time  steps  on  the  acoustic  waves,  while  the 
physics  is  mainly  driven  by  the  solenoidal  part  of  the  flow', 
whose  time  scale,  being  associated  to  advection,  is  larger.  Then 
a three-level  backward  differentiation  formula  is  used  to 
approximate  the  temporal  derivative  of  Qc  leading  to  second- 
order  accuracy.  An  approximate  Newton  method  is  employed  to 
solve  the  non-linear  problem.  At  each  iteration  of  this  inner 
process,  the  inversion  of  the  linear  system  relies  on  Lower- 
Upper  Symmetric  Gauss-Seidel  (LU-SGS)  implicit  method, 
originally  proposed  in  reference  More  details  about  these 
numerical  points  are  available  in  references  n'\ 

Usually  LES  requires  high-order  centered  scheme  for  the 
Euler  fluxes  discretization  (with  spectral-like  resolution)  in 
order  to  minimize  dispersive  and  dissipative  numerical  errors. 
However  such  schemes  cannot  be  applied  easily  in  complex 
geometry.  Indeed,  most  of  aerodynamic  codes  able  to  deal  with 
such  a geometry  are  based  on  finite-volume  technique  in  order 
to  handle  degenerated  cell.  Thus  getting  high-order  method 
becomes  very  time  consuming  due  to  high-order  quadrature 
needed  to  compute  the  fluxes  along  the  cell  boundaries.  As 
several  works  (see  11  for  instance)  have  shown  that  LES  can  be 
carried  out  with  low-order  centered  scheme  in  case  of  sufficient 
mesh  resolution,  only  second-order  accurate  scheme  is 
employed  in  this  study.  But  a special  effort  has  been  carried  out 
to  minimize  the  intrinsic  dissipation  of  the  scheme. 

Thus  the  second  order  accurate  hybrid  upwind/centered 
discretisation,  developped  in  2I,  has  been  used  in  this  study  to 
achieve  a good  compromise  between  robustness  and  accuracy. 
The  key  point  of  the  scheme  is  the  use  of  a sensor,  which  allow's 
to  introduce  locally  some  numerical  dissipation  when  a 
numerical  wiggle  is  detected  on  one  of  the  primitive  variables. 
Therefore,  the  effect  of  the  SGS  model  is  not  affected  by  the 
numerical  dissipation  as  long  as  odd-even  decoupling  is  not 
detected  in  the  flow  solution. 

Geometry,  computational  grid  and  numerical  procedure 

A two-dimensional  (constant  section  in  the  spanwise 
direction)  NACA  0012  airfoil  with  a C=  0.6096  m chord  and  a 
blunted  TE  of  thickness  //  = 2.5  mm  (0,4%  of  the  chord)  has 
been  retained  for  this  study,  with  reference  to  the  airfoil  noise 
experiment  conducted  at  NASA  in  1980  by  Brooks  and 
Hodgson  J\  In  this  experiment,  the  airfoil  had  a span  5 = 0.46  m 
(or  75  % of  the  chord),  whereas  in  the  present  simulation,  the 
LES  computational  domain  has  a spanwise  extent  representing 
only  3.3  % of  the  chord. 


Figure  1 : Computational  grid. 


The  upstream  flow  velocity  is  69.45  m/s  (maximal  velocity 
in  ”),  the  Mach  number  is  0.205  and  the  Reynolds  number  based 
on  upstream  velocity  and  chord  is  2.86  millions.  The  airfoil 
angle-of-attack  with  respect  to  the  flow  direction  is  5°. 

The  3D  curvilinear  computational  grid  is  obtained  by 
replication  in  the  spanwise  direction  (v)  of  a 2D  curvilinear 
structured  grid  made  of  two  domains.  Domain  #1  is  located 
upstream  the  TE  (C-shape),  with  309  points  along  the  airfoil 
body  and  97  points  in  the  radial  direction.  Domain  #2  is  located 
downstream  the  TE,  with  227  points  in  the  (z)  direction 
(including  35  points  on  the  TE  bluntness)  and  103  points  in  the 

(x)  direction.  In  order  to  minimize  the  influence  of  boundary 
conditions,  a large  computational  domain  is  used.  As  shown  in 
Figure  1,  this  2D  grid  extents  nearly  10  chords  above  and  below 
the  profil,  as  well  as  upstream  and  downstream.  Regions  in  the 
vicinity  of  solid  walls  are  highly  refined.  The  whole  mesh  is 
made  of  1 .76  million  points. 

The  horizontal  plane  z = 0 is  a plane  of  symmetry.  In  any 
(x,  z)  plane,  the  smallest  cells  are  located  at  the  TE  comers,  with 
dimensions  (in  wall  units)  Az  = 1 .5  in  the  direction  normal  to 
the  wall  and  Ax'  = 15  along  the  chord.  From  the  comers,  the 
grid  is  streched  in  the  streamwise  direction  (x)  with  a stretching 
coefficient  equal  to  1.06  (upstream)  and  1.09  (downstream).  In 
the  (z)  direction,  the  stretching  coefficients  are  equal  to  1.23 
towards  the  z - 0 plane  and  1.12  towards  the  grid  periphery. 

This  2D  grid  is  replicated  32  times  in  the  spanwise  direction 

(y) ,  with  a constant  step  Ay  = 10  ’C. 

The  mesh  resolution  is  close  to  that  of  Weber  and  Ducros 
who  have  shown  that  such  resolution  is  sufficient  to  obtain 
reasonable  results  for  such  flows. 

A no-slip  condition  is  applied  at  the  airfoil  surface  and  a 
periodic  conditon  is  imposed  in  the  spanwise  direction.  Non- 
reflecting characteristic  boundary  conditions  are  applied  for  the 
far  field.  Moreover,  a steady  RANS  computation  using 
Baldwin-Lomax  models  provides  an  initial  flow  solution.  In 
order  to  ensure  the  time-accuracy  of  the  results,  the  physical 
time  step  is  taken  equal  to  0.5  (is,  meaning  a sampling  frequency 
of  2 MHz.  An  initial  phase  of  100,000  time  steps  was  achieved, 
then  the  useful  computation  was  performed  over  130,000  time 
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steps  with  a total  duration  of  65  ms,  representing  the  convection 
of  the  flow  over  7. 5 airfoil  chords.  The  requirements  of  acoustic 
computations  led  to  a storage  of  one  time  sample  every  100, 
meaning  a storage  sampling  of  20  kHz,  with  a useful  frequency 
band  of  10  kHz.  The  total  computation  cost  was  360  CPU  hours 
on  a NEC-SX5  computer. 


Figure  2 : Contours  of  instantaneous 
Mach  number  isovalues. 


Aerodynamic  results 

Figure  2 displays  contours  of  instantaneous  Mach  number 
isovalues,  showing  the  development  of  the  TBLs  and  the  wake. 
The  transition  between  laminar  and  turbulent  flow's  occurs  at 
x/C  = 0. 1 5 on  the  succion  side  and  at  x/C  =0.9  on  the  pressure 
side.  Mean  velocity  profiles  near  the  TE  were  obtained  by  time- 
averaging instantaneous  data.  The  integration  of  these  velocity 
profiles  provided  TBL  momentum  thicknesses  of  $ = 0.8  mm 
on  the  pressure  side  and  5,  = 1.5  mm  on  the  succion  side. 
Brooks  and  Hodgson's 24  measurement  on  a similar  airfoil  was 
slightly  higher  ( § ~ 3.9  mm),  which  could  be  explained  by  the 
roughness  strips  of  carbon  random  they  used  to  artificially  trip 
the  boundary  layer  at  x/C  = 0. 1 5 . 


Figure  3 : Flow  streamlines  at  the  trailing  edge,  instantaneous 
(above)  and  time-averaged  (below). 

It  is  thougth  that  the  excessive  slenderness  of  the  computed 
TBL  may  be  also  explained  by  an  unsufficient  grid  refinement 


in  the  chord  direction  x near  the  walls  in  the  transition  region. 
However,  it  will  be  shown  in  the  next  part  that  the  computed 
flow  presents  general  qualitative  and  quantitative  features  which 
makes  it  suitable  for  aeroacoustic  applications. 

Figure  3 shows  details  of  the  instantaneous  (above)  and 
time-averaged  (below)  flow  streamlines  around  the  blunted 
trailing  edge.  Time-averaged  results  have  been  favourably 
compared  with  the  results  of  the  mean  flow  computed  separately 
using  a RANS  code. 


UNSTEADY  PRESSURE  FIELD  ANALYSIS 

Surface  pressure  fluctuations 
Time  histories 

Figure  4 (above)  shows  time  histories  of  the  wall  pressure 
fluctuations  near  the  TE  at  symmetrical  points  located  on  the 
pressure  and  succion  sides.  This  plot  exhibits  a time  period  of 
2 ms,  meaning  a dominant  mechanism  at  a frequency  near 
5 kHz,  and  a clear  phase  opposition  between  both  sides.  This 
opposition  is  confirmed  by  the  cross-spectrum  phase  between 
boths  signals  shown  on  Figure  4 (bottom).  These  features  are 
typical  of  a mechanism  of  alternated  vortex  shedding  at  the  TE. 

Power  spectral  densities 

Figure  5 shows  the  evolution  along  the  chord  (from  the 
leading  edge  up  to  the  TE)  of  the  wall  pressure  power  spectral 
density  (0-10  kHz)  on  the  succion  side  (above)  and  the  pressure 
side  (below).  Spectra  are  characterized,  on  both  airfoil  sides,  by 
a narrowband  component  emerging  out  of  a wideband 
continuum. 

Located  near  4.7  kHz,  the  narrowband  bump  is  due  to  the 
vortex  shedding  mechanism.  It  is  perceived  on  the  whole  chord 
length  with  identical  levels  on  both  airfoil  sides.  It  will  be 
shown  later  that  this  component  corresponds  to  acoustic  waves 
generated  at  the  TE  and  propagating  in  the  upstream  direction. 


Figure  4 : Wall  pressure  fluctuations 
on  pressure  and  succion  sides  near  the  TE. 
Above : time  histories.  Below : cross-spectrum  phase. 
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data  at  the  same  position  (x/C=  0.9958).  The  comparison  is 
qualitatively  satisfying  : same  narrowband  component  emerging 
by  =7  dB  from  a wideband  continuum.  However,  the 
quantitative  comparison  shows  that  predicted  levels  are 
overestimated  by  nearly  3 dB,  and  that  the  predicted  vortex 
shedding  frequency  (4.7  kHz)  is  also  overestimated  (3.5  kHz  in 
the  experiment).  Since  this  frequency  is  governed  by  a Strouhal 
number  based  on  the  wake  thickness  (sum  of  the  TE  thickness 
and  TBL  displacement  thicknesses),  this  overestimation  may  be 
partly  explained  by  the  excessive  slenderness  of  the  computed 
turbulent  boundary  layer  (TBL)  with  respect  to  the  experimental 
one.  However  this  cannot  be  solely  explained  by  the  absence  of 
transition  triggering  in  the  simulation. 


Figure  5 : Evolution  of  the  wall  pressure  PSD  along  the  chord. 

The  wideband  continuum  is  generated  by  the  convected 
turbulent  boundary  layers.  Since  the  airfoil  has  a positive  5° 
angle-of-attack,  the  laminar-to-turbulent  transition  occurs  at 
15%  of  the  chord  on  the  succion  side,  and  at  90%  of  the  chord 
on  the  pressure  side.  From  this  transition  point  towards  the  TE. 
the  boundary  layer  thickness  grows,  enriching  spectra  with 
higher  frequency  components.  Notice  that,  on  both  airfoil  sides, 
the  vortex  shedding  component  emerges  from  the  continuum 
near  the  TE  (x/C  > 0.98),  then  is  masked  by  the  nearfield 
boundary  layer,  and  finally  emerges  again  at  x/C  < 0.85. 


Wall  pressure  PSD  at  2.54  mm  from  TE  (x/C  = 0 9958) 


Mwjsuremeal  (toots  ond  Hodgson) 

Compulation  LES  (jycclon  tkn) 

Compulation  L£S  (pressure  side) 

Figure  6 : Wall  pressure  PSD  near  the  TE.  Comparison  of 
simulation  (both  airfoil  sides)  and  experiment 

Comparison  to  experimental  data 

In  their  airfoil  noise  experiment  (performed  with  same 
velocity,  tripped  TBL,  but  zero  incidence)  ”,  Brooks  and 
Hodgson  measured  surface  pressure  with  arrays  of  pinhole 
sensors  of  diameter  0.34  mm.  The  closest  distance  between 
sensors  and  the  TE  was  at  2.54  mm.  Figure  6 compares  the 
spectrum  measured  by  this  sensor  with  simulated  data  computed 


Wavemmher-fretmenrv  spectra 

The  spectral  densities  presented  above  are  single-point  data 
which  do  not  provide  any  information  on  the  propagative 
features  of  the  wall  pressure  field.  This  information  can  be 
revealed  by  processing  a discretized  space-time  Fourier 
transform  of  any  function  p(x,t) , providing  the  so-called 
wavenumber-frequency  spectrum  <t>(k,Co).  For  example,  a 
mono-dimensional  wave  with  frequency  <U0and  velocity  Vu 
propagating  in  the  x direction,  say  p(x,  t)=  e will 

appear  in  the  (k,0))  domain  as  a Delta-function  centered  at 
k = a>JVls  and  ft)  = ft) „ . With  the  same  process,  a white  noise 
propagating  at  celerity  c0  in  the  x direction  will  appear  as  a 
ridge  along  the  k - of  au  line,  or  along  the  k = — 0)1  au  line  if 
the  sound  travels  in  the  opposite  (-x)  direction. 


s moo  ‘woa  booo  tax  inooo 


Frequency  (Hz) 

Figure  7 : Wavenumber-frequency  spectrum  of  wait  presssure 
fluctuations  on  the  airfoil  succion  side  at  x/C  — 0.9. 

Figure  7 presents  the  results  of  this  process  applied  to  the 
wall  pressure  field  computed  via  LES  on  the  airfoil  pressure  side 
near  the  TE  at  x/C=  0.9.  The  wavenumber  domain  extent 
[-tf/A.y,  id  Av]  depends  on  the  local  space  sampling  Ax  of  the 
pressure  field,  or  the  size  of  the  grid  cells.  A wideband 
"convective  ridge"  can  be  observed  along  the  k — of  U t line, 
generated  by  the  convection,  in  the  flow  direction,  of  eddies 
inside  the  TBL  at  an  average  velocity  Uc  slightly  inferior  to  the 
local  mean  flow'  velocity  Uf  . Another  wideband  component  can 


Wavenumber 
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be  observed  along  the  k = -*/(«<>  — Ue  ) corresponding  to  an 
acoustic  propagation  in  the  direction  opposite  to  the  flow,  with  a 
peak  at  4.7  kHz,  the  vortex  shedding  frequency. 
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Figure  8 : Same  as  Figure  7 but  x/C  = 0.5. 


Figure  8 presents  the  same  results  obtained  at  mid  chord 
(x/C  = 0.5).  The  convective  ridge  and  the  upstream  acoustic 
component  are  again  present.  One  interesting  result  is  a new 
component  k = a/(ao  + Ue ) which  corresponds  to  acoustic 
waves  (also  with  a peak  at  4.7  kHz)  propagating  in  the  flow 
direction  (downstream),  which  may  be  generated  by  the 
diffraction,  at  the  airfoil  leading  edge,  of  the  acoustic  waves 
propagating  upstream. 

Flow  pressure  fluctuations 
Instantaneous  maps  of  pressure  fluctuations 

The  existence  of  a sound  field  generated  at  the  TE,  strongly 
suggested  by  the  surface  pressure  results  above,  is  confirmed  on 
Figure  9 which  shows  instantaneous  isovalues  of  the  pressure 
fluctuations  inside  the  flow.  At  every  point  of  the  LES  grid, 
pressure  fluctuations  are  computed  by  substracting  the  time- 
averaged  pressure  to  the  instantaneous  pressure.  The  color  table 
is  adapted  to  the  low  values  of  the  acoustic  fluctuations  (±  6 Pa) 
although  there  are  much  larger  pressure  fluctuations  (up  to 
± 500  Pa)  in  the  TBL  and  in  the  wake.  Concentric  waves  are 
clearly  observed  near  the  TE,  with  a wavelength  corresponding 
to  the  vortex  shedding  frequency  4.7  kHz  which,  as  shown 
before,  dominates  the  source  spectrum  by  a few  decibels.  An 
animation  was  generated  from  300  identical  maps  plotted  for 
successive  time  steps  (A  1=  5.10'!  s)  allowing  to  simultaneously 
visualize,  with  their  respective  celerity,  the  acoustic  radiation  in 
the  fluid  and  the  turbulence  convection  along  the  airfoil  chord. 


Figure  9 : Instantaneous  isovalues  of  pressure 
fluctuations  obtained  from  LES  data. 


Figure  10  : Evolution  of  the  wall  pressure  spectra  along  the 
vertical  grid  line  x - C (starting  from  the  TE  upper  comer)  with 
respect  to  the  vertical  distance  z. 


Low-nass  filtering  bv  the  CFD  grid 

It  is  interesting  to  notice  that  the  wave  pattern  corresponding 
to  the  vortex  shedding  noise  vanishes  at  a half-chord  from  the 
airfoil,  when,  at  the  same  time,  larger  wavelength  are  observed 
much  farther.  This  is  explained  by  the  radial  stretching  of  the 
grid  LES  which  acts  on  the  noise  field  as  a low-pass  filter.  It  is- 
known  that  the  propagation  of  an  acoustic  wave  will  not  be 
correctly  simulated  if  it  is  discretized  by  using  less  than  4 or  6 
cells  by  wavelength. 

To  illustrate  this.  Figure  10  shows  the  evolution  of  the  wall 
pressure  spectra  along  the  vertical  grid  line  x = C (starting  from 
the  TE  upper  comer)  with  respect  to  the  distance  from  the  TE 
(in  logarithmic  scale).  This  plot  clearly  shows  that  low 
frequency  waves  propagate  much  farther  than  high  frequency 
waves. 

Figure  1 1 shows  the  same  information  for  several  frequency 
bandwidths  in  which  spectra  are  integrated.  Each  frequency 
bandwidth  shows  similar  behaviour  : (i)  at  short  distance  from 
the  TE,  spectra  are  dominated  by  high-level  turbulence  ; (ii) 
then  there  exists  a domain  where  the  pressure  fluctuations 
decrease  as  s'  , which  is  typical  of  a two-dimensional  acoustic 
field  ; (iii)  beyond  a cut-off  distance  the  levels  vanish.  It  was 
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checked  that,  for  any  frequency  bandwidth,  this  cut-off  distance 
corresponds  to  a grid  location  where  cells  become  larger  that  a 
quarter  acoustic  wavelength. 


Figure  11 : Same  information  as  in  Figure  10, 
presented  in  integrated  frequency  bandwidth s 
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Figure  12  : Same  information  as  in  Figures  10  and  11,  but 
presented  through  a wavenumber-frequency  spectrum  . 


Figure  12  shows  .a  wavenumber  frequency  spectrum 
computed  from  the  LES  pressure  fluctuations  field  along  the 
same  vertical  line  x = C.  The  ridge  located  along  the  "acoustic 
line"  k = ft/  a0  clearly  demonstrates  the  radiative  acoustic 
nature  of  this  pressure  field  with,  again,  a peak  of  power  at  the 
vortex  shedding  frequency  4.7  kHz.  A much  weaker  component 
is  also  visible  on  the  line  k = correponding  to  waves 

propagating  towards  the  TE,  but  this  peak  lays  more  than  25  dB 
below'  the  radiative  component  and  may  even  not  be  physical 
with  regards  to  the  resolution  of  the  process. 


Snanwise  coherences 

Inside  the  TBL  and  the  wake,  turbulent  structures  or  eddies 
are  convected  at  low-  velocities  (<  (/,,)  which  means  that  typical 
spanw'ise  lengh  scales  may  be  smaller  than  the  spanwise  extent 
of  the  LES  domain,  and  correctly  simulated.  To  illustrate  this, 
Figure  13  shows  the  spanwise  evolution  of  the  coherence 


y2(Ay,/)  of  the  surface  pressure  field  near  the  TE  in  several 
frequency  bandwidths.  y ~ (Ay,  /)  is  defined  as  : 


72(  Ay,/)  = 


M0./)|s„(Ay,/)| 


m 


where  Sp(v,f ) denotes  the  surface  pressure  spectrum  at  the 
spanwise  position  y and  S (yl.y2,f)  denotes  the  surface 

pressure  cross  spectrum  between  spanwise  positions  v,  and  v,. 
The  coherence  rapidly  vanishes  at  any  frequency  bandwidth, 
meaning  that  integral  length  scales  are  smaller  than  the  spanwise 
extent  of  the  LES  domain.  It  can  be  seen  that  the  largest  length 
scales  correspond  to  the  vortex  shedding  frequency,  a 
mechanism  which  significantly  correlated  in  the  spanwise 
direction.  It  should  be  noted  that  these  spanw'ise  coherences 
provide  wall  pressure  spanwise  integral  length  scales  which  are 
favorably  compared  to  experimental  values  measured  by  Brooks 
and  Hodgson  on  their  NACA0012  airfoil  w. 


0.2  - 2.0  kit  - Center  1.1  kHz 
2.0  - S£  kHz  - Center  2,9  kHz 
3.a  - 5.6  kHz  - Center  4.7  kHz 
5.6  - 7.4  kHz  - Ce 6 5 kHr 
7.4  - 0.2  kHr  - Center  0 3 kHz 


Figure  13  : Spanwise  evolution  of  the  coherence  of  the  surface 
pressure  field  on  the  suction  side  near  the  TE  (x/C  = 0.9958). 
Integrated  frequency  bandwidths. 


Figure  14  show's  spanwise  coherence  of  the  pressure  field 
above  the  TE  on  the  grid  line  j„  = 54  located  at  distance  zu  = 37.9 
mm  from  the  TE,  where  pressure  fluctuations  are  expected  to  be 
mostly  acoustic  : as  explained  below,  this  particular  grid  line 
will  be  chosen  for  the  generation  of  the  integration  surface  used 
in  the  KirchhofT  prediction.  Figure  14  shows  that  typical  length 
scales  seem  to  be  much  longer  than  the  spanwise  extent  of  the 
LES  grid,  especially  at  low  frequency  and  the  frequency  of 
vortex  shedding. 

In  other  words,  the  very  nearfield  turbulent  flow  is  expected 
to  be  strongly  three-dimensional,  whereas  the  acoustic  field 
simulated  at  some  distance  of  the  airfoil  looks  two-dimensional. 
This  result  will  have  significant  consequences  in  the  process  of 
spanwise  replication  of  LES  data. 
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Figure  14  : Spanwise  evolution  of  the  coherence  of  the  pressure 
field  above  the  TE  on  the  grid  line  j0  = 54  at  x/C=  1.  Same 
integrated  frequency  bandwidth:  as  Figure  14. 


KIRCHHOFF  ACOUSTIC  PREDICTION 

Introduction 

The  preceding  results  show  that  the  compressible  LES  has 
the  capacity  of  simulating  the  generation  of  the  aeroacoustic 
mechanism  of  TE  noise,  including  both  aerodynamic  sources 
and  local  radiated  field.  Unfortunately,  due  the  strong  grid 
stretching,  this  acoustic  field  rapidly  collapses  at  rather  short 
distance  from  the  airfoil.  As  a consequence  of  the  last  part,  it 
makes  sense  to  find  a surface,  enclosing  the  airfoil,  close  enough 
to  the  body  to  ensure  that  the  simulated  acoustic  field  has  not 
been  low-pass  filtered  by  the  grid,  or  at  least  still  contains  the 
low-frequency  bandwidth  that  we  are  interested  in.  Moreover,  if 
this  surface  can  be  proved  free  of  any  vorticity  convected  by  the 
flow,  it  is  a perfect  candidate  for  the  application  of  a Kirchhoff 
integration. 

Formulation 

A 3D  Kirchhoff  formulation  is  implemented  in  the  KIM 
(Kirchhoff  Integration  Method)  code  which  was  developed  at 
ONERA  for  helicopter  rotor  noise  prediction  ” and  jet  noise 
prediction  2\  This  code  provides  the  noise  radiated  by  any  three- 
dimensional  surface  in  a flow  with  uniform  velocity  (/„  in  the 

direction  X1 , given  the  density  (or  pressure,  assuming  isentropy 
2 

or  p'  = dip' ) and  its  normal  gradient  along  this  surface. 


The  formulation  is 21 : 


(10) 


in  which  Dxr  is  the  distance,  corrected  from  convection  effects, 
between  the  source  point  ¥ , located  on  the  integration  surface, 
and  the  observer  point  X (both  defined  in  a reference 
framework  in  which  the  surface  S is  at  rest) : 


®xr  — 


d-M^-Y,) 

P2 


(U) 


fK{x,  P,  r) is  given  by  : 


fK(x . ?.*)  = ■$■(!- <>■(* - P)  + <«,  §7 - ^ 


A/0n,+ 


dr 


(12) 


This  formulation  assumes  that  the  pressure  field  on  S 
satisfies  the  convected  wave  equation  : 


(13) 

which  governs  the  propagation  of  acoustic  waves  in  a medium 
with  uniform  velocity  t/„  in  the  direction  X, . 


Snace-time  discretization 

The  space-time  discretization  of  the  Kirchhoff  integration  on 
a given  closed  control  surface  assumes  that  temporal 
fluctuations  of  fluid  pressure  are  time-sampled  (with  time  step 
At  and  duration  N&t)  at  any  cell  (with  maximal  dimension  Act) 
of  a 2-layer  (with  maximal  normal  separation  An  between  both 
layers)  surface  grid.  If  A,  and  f = XJc0  denote  the  acoustic 
wavelength  and  period,  the  crucial  parameters  of  the  process  are 
Act  /Xa,  An/ Xa  and  At/T.  These  parameters  were  carefully 
adjusted  on  the  basis  of  a parametric  study  (that  can  be  found  in 
the  Annex  of  Ref.  **)  which  was  performed  before  to  achieve 
actual  Kirchhoff  integration. 


Construction  of  the  control  surface 

Figure  15  shows  how  the  integration  surface  is  constructed 
from  one  given  C-shaped  grid  line  surrounding  the  airfoil  in 
domain  #1,  extended  in  domain  #2  by  the  two  corresponding 
horizontal  grid  lines.  These  horizontal  lines  are  stopped  at  nearly 
one  third-chord  from  the  TE,  where  acoustic  levels  are  very  low 
and  grid  cells  become  much  larger  than  the  acoustic  wavelength. 


Figure  15  : Kirchhoff  integration  : 
construction  of  the  control  surface. 


Then  the  three-dimensional  surface  is  obtained  by 
N replications  in  the  spanwise  direction  along  the  desired  span 
extent  L = NAL.  In  the  present  application,  with  reference  to 
Brooks  and  Hodgson's  experiment 2,!,  the  airfoil's  span  was  taken 
equal  to  L = 0.46  m or  0.75C.  N was  fixed  to  N = 92  according 
to  a parametric  study  presented  in  the  Annex  of  Ref. 2". 

It  should  be  noted  that  the  final  3D  surface  is  not  closed, 
neither  at  its  two  lateral  sides  (their  contribution  to  the  noise 
radiated  in  the  z direction  is  neglectible)  nor  at  its  downstream 
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extremity  (which  contains  very  low  acoustic  levels  and  may  be 
crossed  by  undamped  vorticity  convected  in  the  wake). 

The  LES  3D  grid  has  a spanwise  extent  limited  to  3.3  % of 
the  airfoil  chord,  so  the  LES  data  have  to  be  extrapolated  in  the 
spanwise  direction.  The  2D  behaviour  of  the  simulated  acoustic 
field,  especially  confirmed  by  the  spanwise  coherences 
presented  on  Figure  14,  suggest  to  simply  replicate,  along  the 
whole  span  L,  the  2D  field  computed  at  y = 0. 

Finally  the  critical  design  parameter  of  the  integration 
surface  is  the  index  jn  of  the  grid  line  surrounding  the  airfoil, 
which  varies  from  j0  = 1 on  the  airfoil,  to  ju  - 97  at  the  grid 
border.  Once  this  parameter  is  chosen,  the  normal  separation  An 
is  directly  related  to  j,„  as  grid  lines  of  indexes  j0  and  j„  +1  are 
used  to  compute  the  pressure  normal  gradients. 


Figure  16  : Parametric  study  of  the  position  of  the  Kirchhoff 
integration  surface  S (parameters  : index  ja  of  the  used  grid  line, 
and  normal  distance  z„  to  the  airfoil  TE). 
Noise  radiated  atz0  = 2C  above  the  trailing  edge. 


Figure  16  shows  the  influence  of  jt  on  the  power  spectral 
density  of  the  noise  radiated  at  z = 2C  above  the  trailing  edge, 
computed  by  the  Kirchhoff  integration.  It  should  be  noted  that 
all  other  parameters  were  adjusted  to  final  values  deduced  from 
the  parametric  studies  described  in  the  Annex. 

On  each  presented  spectrum,  the  index  y'0and  the  distance  of 
the  grid  line  from  the  airfoil  are  mentioned.  This  plot  shows  that 
a family  of  intermediate  grid  lines  provide  very  similar  spectra, 
which  is  a good  indication  of  a region  where  the  pressure  field  is 
mostly  acoustic.  The  spectra  computed  from  the  grid  lines 
nearest  to  the  airfoil  are  obviously  contaminated  by  intense 
vorticity  convected  in  the  TBL  and  the  wake.  On  the  other  hand, 
spectra  computed  from  farther  grid  lines  are  progressively 
attenuated  in  the  high  frequency  band,  which  is  another 
confirmation  of  the  low-pass  filtering  generated  by  the  grid 
stretching.  The  grid  line  ja  = 54  located  at  the  distance 
z„  = 37.9  mm  from  the  airfoil  was  finally  chosen  for  the 
Kirchhoff  computations.  The  normal  distance  between  lines  54 
and  55  is  An  = 4.7  mm. 


noise  can  be  computed  at  thousands  of  locations  to  build 
instantaneous  maps  or  time-animated  pressure  distributions. 


Figure  J 7 ; Noise  contours  computed  via  the  Kirchhoff  method 
outside  the  integration  surface. 


Figure  18:  Pressure  contours  inside  (LES)  and  outside 
(Kirchhoff)  the  integration  surface 


Figure  17  shows  such  a map  computed  on  a two- 
dimensional  cartesian  grid  of  regular  space  step  Ax  = Ay  = 1.2 
cm  which  extends  two  chords  around  the  airfoil.  It  is  interesting 
to  compare  this  Kirchhoff  noise  field  to  the  similar  result  of 
Figure  9,  which  was  obtained  from  the  LES  pressure  fluctuation 
field.  Both  figures  are  using  the  same  colour  scale,  showing  how 
the  acoustic  field  is  rapidly  filtered  and  eliminated  from  the  LES 
data. 

Figure  18  provides  a closer  view  of  the  control  surface 
region,  in  order  to  evaluate  the  continuity  between  the  pressure 
fields  inside  and  outside  the  surface.  Inside  the  surface,  the 
pressure  field  is  directly  deduced  from  the  LES  data  as  in  Figure 
9.  Outside  the  surface,  the  pressure  is  computed  via  a Kirchhoff 
integration  as  in  Figure  17. 


Results 

Instantaneous  pressure  distributions 

The  Kirchhoff  integration  method  offers  the  conveniency  of 
easily  computing  the  noise  at  any  location  outside  the  control 
surface.  If  the  available  computing  resources  are  sufficient,  the 


Farfield  noise 

Figure  16  shows  noise  spectra  computed  from  various 
Kirchhoff  surfaces  (including  the  surface  based  on  the  grid  line 
j0  = 54  z0  = 37.9  mm  which  was  finally  chosen)  at  z = 2 C.  This 
specific  position  was  chosen  because  it  corresponds  to  the 
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location  of  a microphone  in  the  Brooks  and  Hodgson's 
experiment  :4.  The  PSD  measured  by  this  microphone  is 
presented  on  Figure  20  (dashed  line).  The  noise  predicted  via 
the  Kirchhoff  integration  overestimates  the  measured  data  by 
more  than  15  dB.  This  overestimation  may  have  several  reasons 
acting  together.  Firstly,  it  is  important  to  distinguish  the 
nearfield  unsteady  CFD  and  the  farfield  noise  prediction  : it  has 
been  seen  that  wall  pressure  fluctuations  already  overestimated 
measured  values  by  3 to  5 dB,  which  means  that  the 
overprediction  by  the  acoustic  prediction  may  be  closer  to  10 
dB. 

Secondly,  the  spanwise  replication  of  the  central  pressure 
field  assumes  that  the  acoustic  field  on  the  control  surface  is 
perfectly  correlated  (or  fully  two-dimensional)  in  the  spanwise 
direction.  This  assumption  was  considered  in  accordance  with 
the  spanwise  coherences  shown  on  Figure  14.  However,  it  is 
almost  certain  that  these  coherences  are  artificially 

overestimated  due  to  the  small  extent  of  the  LES  computational 
domain  (in  a similar  approach  based  on  the  Ffowcs  Williams 
and  Hall's  TE  noise  theory  applied  to  incompressible  LES, 
Wang  * showed  that  the  radiated  noise  computed  from  CFD  data 
replicated  in  the  spanwise  direction  is  very  sensitive  to  the 
spanwise  coherence  of  these  data). 


FFOWCS  WILLIAMS-HA WRINGS  INTEGRATION 


Formulation 

Since  wall  pressure  spanwise  coherences  (Figure  13)  and 
fluctuation  levels  (Figure  6)  were  found  to  match  experimental 
data,  a second  noise  prediction  was  achieved  using  the  Ffowcs 
Williams-Hawkings  (FW-H)  method  (dipole  surface  term  only). 

A 3D  FW-H  formulation  is  implemented  in  the  FIM  (FW-H 
Integration  Method)  code  which  was  developed  at  ONERA  for 
helicopter  rotor  noise  prediction.  This  code  provides  the  noise 
radiated  by  any  three-dimensional  surface  in  a flow  with 
uniform  velocity  Ua  in  the  direction  Xt , given  the  data  p,  p and 
puj  along  this  surface.  The  formulation  isw  : 

= + j dSdr  (14) 


where  Dxr  and  d are  defined  in  ( 1 1 ) and  : 


J d_ 

p2  dr 


l(x-F)h 


and 


+±[p2a.(x-y)+(m^  -(y0B)(x,  -k,)] 

(15) 
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and  all  other  notations  are  identical  as  in  the  Kirchhoff 
formulation.  It  should  be  noted  that  this  formulation  is  not 
complete  since  it  ignores  the  volume  integration  term  : when  the 
surface  is  taken  at  the  airfoil  wall,  the  volume  outside  the 
surface  obviously  contains  a distribution  of  quadrupoles,  but  it  is 
assumed  that  this  term  can  be  neglected  (this  assumption  has 
been  verified  in  aeroacoustic  mechanisms  involving  solid  walls, 
such  as  helicopter  rotor  noise  or  propeller  noise). 

Application 

The  application  of  the  above  formulation  to  the  present  case 
relies  on  the  comparison  of : 

• the  spanwise  integral  length  scale  of  wall  pressure 
fluctuations  computed  by  LES  : Ay  <3  mm, 

• the  spanwise  extent  of  the  LES  computational  domain : 
Sia  ~ 20  mm  or  3,3  % of  the  chord, 

• the  spanwise  extent  of  the  real  airfoil  used  in  the 
experiment  S.:XF  = 0,46  m. 

In  a first  step,  the  FIM  code  is  used  to  compute  the  DSP 
0Lt:/(o)  of  the  noise  radiated  by  a narrow  airfoil  of  span  5,tr 
This  is  done  by  using  a surface  located  at  the  airfoil  wall  and 
integrating  the  data  computed  via  LES  on  all  the  33  planes  of 
the  3D  LES  domain.  It  is  then  assumed  that  the  real  airfoil  of 
span  S,„,  is  made  of  SafSus  narrow  airfoils  located  side-by-side. 
Since  Ar  < St£r  any  couple  of  these  narrow  airfoils  behave  as 
fully  decorrelated  sources  and  an  estimation  of  the  DSP  &FJI/(0) 
of  the  noise  radiated  by  the  real  airfoil  of  span  Sc„  is  given  by  : 

®EXp(<°)=^XL't>LEs(a>)  (17) 

‘bLES 

Results 

Figure  20  shows  the  noise  spectrum  estimated  via  this 
numerical  procedure  at  z = 2C  above  the  TE.  It  is  compared, 
without  correction  or  normalisation,  to  the  experimental  data 
obtained  by  Brooks  and  Hodgson  at  the  same  position. 

It  is  interesting  to  notice  that  experiment  and  simulation 
differ  in  the  same  way  (frequency  shift  and  amplitude 
difference)  as  wall  pressure  spectra  shown  on  Figure  6,  which 
may  confirm  the  validity  of  the  acoustic  method. 


Noim  spectro  nt  R = ) 2 r or  2C 


Figure  20  : Power  spectral  density  of  the  noise  radiated  at 
z = 2C  above  the  TE.  Comparaison  of  computed  (solid  line) 
and  measured  (dashed  line)  data  24. 


Directivity  diagrams 

Figure  21  shows  a diagram  obtained  by  computing  (with 
both  Kirchhoff  and  FW-H  methods)  the  noise  radiated  on 
equally  spaced  points  located  on  a circle  of  radius  r=  1.2  m (or 
2C)  centered  on  the  airfoil  TE.  The  directivity  is  compared  to 
the  theoretical  directivity  of  the  noise  emitted  by  turbulence 
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passing  the  trailing  edge  of  an  infinite  half-plane.  It  is  also 
compared  to  experimental  values  given  by  Brooks  and 
Hodgson 


Figure  21  : Directivity  diagram  of  radiated  noise  on  a circle  o f 
radius  2C  centered  on  the  TE.  Comparison  of  computed  data 
(Kirchhoff  and  Ffowcs  Williams  Hawkings  methods)  with  half- 
plane TE  noise  and  experimental  data  * 


LINEARIZED  EULER  EQUATIONS 

The  need 

Going  back  to  the  pressure  contours  shown  on  Figure  18 
obtained  via  the  Kirchhoff  integration,  it  is  interesting  to  follow 
the  acoustic  waves  emitted  at  the  TE  and  propagating  in  the 
upstream  direction,  especially  along  the  pressure  side  chord  on 
which  the  TBL  is  only  visible  near  the  TE.  At  the  vicinity  of  the 
TE,  these  waves  are  rather  continuous  in  both  LES  and 
Kirchhoff  pressure  fields.  When  they  have  traveled  along  more 
than  a half  chord,  it  becomes  clear  that  waves  in  the  LES 
domain  (inside  the  surface)  are  slightly  slower  than  their 
counterpart  propagating  in  the  Kirchhoff  domain  outside  the 
surface.  To  explain  this,  one  must  remember  that  the  Kirchhoff 
method  assumes  that  the  flow  outside  the  controle  surface  has  a 
uniform  velocity  U0,  whereas  the  mean  flow  simulated  via  LES 
inside  the  control  surface  present  local  velocities  Uf  that  may 
be  significantly  larger  than  U0,  especially  on  the  succion  side  of 
a lifting  airfoil.  This  last  result  underlines  what  is  probably  the 
most  critical  limitation  of  integral  methods  such  as  the 
Kirchhoff  integration  : outside  a given  surface  or  volume,  they 
do  not  take  into  account  the  inhomogeneities  of  the  mean  flow. 

This  limitation  fully  justifies  the  development  of  a 
complementary  tool  based  on  the  discretization  of  LEE,  which 
have  the  ability  to  account  for  local  mean  velocity  gradients. 
However,  since  the  application  of  such  methods  may  become 
very  computationally  consuming,  especially  in  complex  3D 
cases,  it  should  be  strictly  limited  to  the  regions  where  strong 
velocity  gradients  make  it  really  necessary.  In  the  present 
problem,  such  region  could  be  reasonably  limited  to  a domain 
extending  one  chord  away  from  the  airfoil.  Outside  this  new 
control  surface,  the  mean  flow  can  be  considered  as  rather 
uniform  and  the  Kirchhoff  method  can  be  applied  to  provide 
accurate  prediction  of  the  sound  at  much  larger  distance  from 
the  airfoil. 


The  E3P  code 

The  E3P  code  (Propagation  via  Euler  Equations  under  a 
Small  Perturbation  Hypothesis)  based  on  discretized  LEE  is 
currently  developped  at  ONERA  for  this  purpose IJ.  The 
formulation  is  based  on  the  splitting  of  the  total  field  in  a mean 
flow  and  a perturbation  field.  The  equations  are  discretized  in 
conservative  form  and  can  be  either  linearized  (only  first  terms 
are  kept)  or  not  (all  linear  and  non-linear  terms  are  used). 
Cartesian  or  curvilinear  2D/3D  (monodomain)  structured  grids 
can  be  used.  In  the  space  domain,  the  discretization  uses  high 
order  finite  differente  explicit  or  implicit  schemes.  In  the  time 
domain,  multi-steps  schemes  (Adams-Bashford  or  Runge-Kutta) 
are  used.  Specific  boundary  conditions  are  implemented  for 
solid  surfaces,  non-reflexive  border,  periodicity  or  data 
injection.  Non-uniform  mean  flows  are  taken  into  account  The 
next  part  shows  several  configurations  of  acoustic  scattering 
problems  solved  with  this  code.  The  last  paragraph  explains  how 
this  code  could  be  coupled  with  unsteady  CFD  computations. 

Acoustic  scattering  using  E3P 

Figure  22  shows  the  free-field  radiation  of  an  acoustic 
harmonic  Gaussian  source  in  a steady  uniform  flow  at  Mach  0.5, 
computed  on  a 200  x 200  weavy  grid  and  compared  with  the 
analytic  solution. 


Figure  22  :Radiation  of  a harmonic  Gaussian  source 
in  a uniform  flow  M = 0.5.  E3P,  weavy  grid 
200  x 200,  compared  with  analytic  solution. 


Figure  23  shows  the  acoustic  scattering  of  a harmonic 
source  on  a rigid  circular  cylinder  of  radius  a in  a fluid  at  rest. 
The  RMS  pressure  field  (amplitude  and  phase)  along  a circle  of 
radius  2a  is  compared  with  results  obtained  with  a Boundary 
Element  Method  (BEM).  Figure  24  shows  the  same 
configuration,  but  with  a non  uniform  (potential)  mean  flow 
around  the  cylinder.  Figure  25-shows  the  acoustic  scattering  of  a 
harmonic  source  on  a rigid  elliptic  cylinder. 

Figure  26  shows  the  acoustic  scattering  of  a harmonic 
source  on  a rigid  Joukowski  airfoil  in  a potential  mean  flow. 
Figure  27  displays  three-dimensional  computation  of  the 
acoustic  scattering  of  a harmonic  source  on  a 3D  rigid  circular 
cylinder,  again  and  with  a potential  mean  flow. 


Figure  24  : Acoustic  scattering  of  an  harmonic  source  on  a 
rigid  circular  cylinder  with  potential  non-uniform  mean  flow 
(. M = 0.3).  RMS  pressure  field. 


The  coupling  of  local  unsteady  CFD  with  discretized  LEE  is 
a complex  problem.  Up  to  now,  most  applications  were  based 
on  die  formulation  of  volume  distribution  of  aeroacoustic 
sources  constmcted  from  the  CFD  results  and  used  as  right- 
hand-side  source  terms  in  Euler  equations. 

In  the  present  application,  it  was  shown  that  compressible 
LES  is  able  to  simulate  the  generation  and  local  propagation  of 
noise.  This  result  suggests  that  LES  Gould  be  coupled  with  the 
noise  E3P  code  using  a surface  interface  similar  to  the  Kirchhoff 
control  surface.  This  basic  principle  requires  that  a perturbation 
field  can  be  injected  at  one  border  of  a LEE  domain,  a 
functionality  which,  in  a first  step,  was  tested  by  injecting  an 
acoustic  field  radiated  by  an  harmonic  source  located  outside  the 
LEE  domain. 

Figure  28  shows  such  validation  in  the  case  of  a a 
cylindrical  grid.  The  (virtual)  monopole  is  located  inside  the 
inner  border  of  the  domain,  and  its  (analytical)  field  is  injected 
on  this  inner  border.  The  acoustic  field  computed  in  the  LEE 
domain  is  favorably  compared  with  the  analytical  field. 


Figure  25  : Acoustic  scattering  of  an  harmonic  source  on  a 
rigid  elliptic  cylinder  (M  = 0).  Comparison  of  amplitude  RMS 
pressure  field  (LHS))  with  BEM  results  (RHS). 
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Figure  28  : Injection,  at  the  inner  border  of  a cylindrical  grid, 
of  the  acoustic  field  radiated  by  a (virtual)  monopole  located 
inside  the  inner  border.  Comparison  with  the  analytical  field. 
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On  Figure  29  the  same  analytical  acoustic  field  is  injected  at 
the  inner  border  of  a grid  which  was  derived  from  the  LES  grid 
shown  on  Figure  1.  Again,  the  pressure  field  computed  in  the 
LEE  domain  is  favourably  compared  with  the  abanlytical  direct 
propagation.  The  next  step  will  consist  in  injecting  LES  data  at 
the  same  inner  border. 

It  is  believed  that  the  next  significant  results  of  CAA  applied 
to  airfoil  noise  predictions  will  be  obtained  via  the  combination 
of  (i)  unsteady  CFD  (LES)  for  the  nearfield  flow,  (ii)  the 
linearized  Euler  equations  for  intermediate  acoustic  sound  field 
and  (iii)  an  integral  method  (Kirchhoff,  Ffowcs  Williams- 
Hawkings)  for  the  prediction  of  the  far  field  noise. 


Figure  29  . Injection,  at  the  inner  border  of  a CFD  grid,  of  the 
acoustic  field  radiated  by  a (virtual)  monopole  located  inside 
this  inner  border.  Comparison  with  the  analytical  field. 


CONCLUSIONS 

The  unsteady  flow  past  an  isolated  NACA0012  airfoil  with  a 
constant  section  and  a blunted  TE  is  computed  via  compressible 
LES.  Detailed  analysis  of  the  results  shows  that  the  unsteady 
flow  is  accompanied  by  an  acoustic  field  generated  by  the 
acoustic  scattering  of  (i)  the  turbulent  boundary  layers 
convected  on  both  airfoil  sides  (broadband  noise)  and  (ii)  the 
alternated  vortex  shedding  generated  by  the  TE  bluntness 
(narrow  band  component).  A close  examination  of  surface 
pressure  fluctuations  near  the  TE  shows  that  spanwise  integral 
length  scales  are  probably  smaller  than  the  computational 
domain  extent  in  that  direction,  although  the  latter  only 
represents  3.3  % of  the  chord.  Surface  pressure  spectra  are 
compared  with  published  experimental  data  obtained  on  the 
same  airfoil  geometry.  The  computed  vortex  shedding  frequency 
is  slightly  higher  than  the  measured  one,  which  can  be  explained 
by  an  underestimation  of  the  TBL  thicknesses.  Predicted  levels 
also  slightly  overestimate  measured  levels.  Due  to  the  strong 
stretching  of  the  LES  computational  grid,  which  acts  as  an 
acoustic  low-pass  frequency  filter,  this  acoustic  field  cannot 
radiate  farther  than  a half-chord  from  the  body. 

In  a second  step,  the  LES  has  been  relayed  by  an  acoustic 
propagation  method  to  simulate  the  farfield  noise.  For  this 
purpose,  two  different  integral  methods  have  been  used,  a 
Kirchhoff  integration  and  a Ffowcs  Williams-Hawkings 
integration. 

The  Kirchhoff  control  surface  is  constructed  via  the 
spanwise  development  of  a grid  line  surrounding  the  airfoil. 
This  method  necessitates  a careful  study  of  several  critical 
parameters  such  as  the  position  of  the  control  integration  surface 
enclosing  the  airfoil  and  its  wake,  the  spanwise  replication  step, 
the  distance  between  the  two  layers  of  the  surface  (for  correct 
computation  of  pressure  normal  gradients)  and  the  time  step,  all 
parameters  being  normalized  by  the  acoustic  wave  length  or 


time  period.  Since  the  acoustic  field  generated  at  the  TE  shows 
strong  two-dimensional  characteristics,  the  pressure  field  on  the 
integration  surface  is  simply  obtained  by  replicating  the  LES 
data  in  the  spanwise  direction.  Predicted  farfield  noise  spectra 
are  compared  with  published  experimental  data  obtained  in  an 
anechoic  facility  with  the  same  airfoil  geometry  (but  much 
larger  span).  Predicted  levels  overestimate  experimental  spectra, 
which  is  most  probably  due  to  the  assumption  that  the  acoustic 
field  is  fully  correlated  along  the  spanwise  extent  of  the 
integration  surface,  which  is  obviously  only  true  at  very  low 
frequencies  and  in  the  frequency  domain  of  the  vortex  shedding. 

In  a second  step,  the  Ffowcs  Williams-Hawkings  method 
was  used  along  with  wall  pressure  data  which  provided 
spanwise  correlations  matching  experimental  values.  In  this 
case,  the  far  field  noise  spectrum  is  obtained  in  two  steps.  In  a 
first  step,  the  noise  radiated  by  an  narrow  airfoil  is  obtained  by 
integrating  wall  pressure  data  in  the  whole  LES  computation 
domain.  Then  the  noise  radiated  by  an  airfoil  with  realistic  span 
is  obtained  by  summing  decorrelated  contributions  of  several 
identical  narrow  airfoils  lying  side-by-side.  Final  results 
favorably  agree  with  experimental  measurement. 

' Finally,  close  examination  of  acoustic  waves  propagating  (i) 
in  the  LES  pressure  field  inside  the  controle  surface  and  (ii)  in 
the  Kirchhoff  pressure  field  outside  the  surface,  clearly  shows 
the  influence  of  the  uniform-flow  assumption  which  is  implicit 
in  most  integral  methods.  This  strongly  justifies  the  use  of  a 
novel  acoustic  method  based  on  the  discretized  LEE.  This 
method  is  currently  under  development  Several  solutions  of 
acoustic  scattering  problems  were  solved  with  this  code  and  are 
presented.  Then  it  is  shown  how  this  tool  could  be  combined 
with  nearfield  LES  on  the  basis  of  a 2D  interface.  It  is  believed 
that  this  combination  will  soon  bring  significant  progress  in  the 
field  of  the  numerical  prediction  of  airframe  noise. 
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